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Single-charge tunneling in ultrasmall voltage biased SIS junctions in a high-impedance elec- 
tromagnetic environment is considered. The Cooper pair current is calculated at T — on the 
basis of the elementary tunnel Hamiltonian for quasiparticles. Therefore, the transfer of Cooper 
pairs emerges automatically as an effect of higher order perturbation theory. The supercurrent also 
depends on the dissipative part of the Josephson current amplitude. 
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I. INTRODUCTION 



Effects of single-charge tunneling in ultrasmall capacitance junctions has become of much experimental and theo- 
retical interest. For a review see for instance Ref. |l|. Modern nanolithography allows the fabrication of junctions with 
capacitances C < 10 _16 F where the electrostatic energy differences dominate thermal fluctuations at the IK scale. 
This opens the door to a new kind of electronics. 

The current through a SIS junction (see Fig. [I]) can be carried by quasiparticles (quasiparticle current (I) qp ) and 
by Cooper pairs (supercurrent (I) s )- Here we develop a nonlinear response approach to the Cooper pair tunneling in 
voltage biased ultrasmall SIS junctions basing on the elementary tunneling Hamiltonian for quasiparticles. Since this 
Hamiltonian describes only the tunneling of single quasiparticles (le) the Cooper pair tunneling (2e) corresponds in 
perturbatkm theory to a process of higher order. In leading order this approach yields the known result of quasiparticle 
tunnelingB. Because of the voltage biasing and due to Eq. (||) the ordinary Josephson current which is also an effect 
of first order perturbation theory does not play any role. 

The current-voltage characteristic of small tunnel junctions is essentially influenced by the external circuit. This 
electromagnetic environment is able to absorb energy which is for superconducting electrodes at zero temperature the 
only possibility to transfer the energy gain of the tunneling process. Because Cooper pairs live in the condensate they 
cannot absorb this energy. To simplify matters we restrict ourselves to the limiting cases of low- and high-resistance 
environments (Re <€i Rq and Re ^> Rq) at T = where the quasiparticle currents are suppressed for voltages lower 
than the thresholds 2A/e and (2A + E c )/e respectively. Rq = h/e 2 is the quantum resistance and 2A labels the 
superconducting energy gap. The additional part E c — hu> c — e 2 /(2C) corresponds to the Coulomb energy. For single 
junctions is known that the Coulomb blockade can only be observed if the junction is sufficiently decoupled from the 
voltage bias by a high-resistance environment. Beyond the thresholds Cooper pairs can break up into quasiparticles 
and the tunnel current is carried mainly by quasiparticles. p, m 

Cooper pair tunneling is described in literature by using the model of an effective HamiltonianEnD with the pertur- 
bation term 

H T = Ejcos& (1) 

where the operator exp(±itf') changes the macroscopic charge Q on the junction by the value ±2e corresponding to 
the charge of a Cooper pair. This means that simultaneously tunneling of two electrons (Cooper pair) is introduced 
by hand. Then the calculated supercurrent readsEm 

neEj f (2eV\ . ( 2eV\ 



with 



P'M = 77- / e 4J(r)+ ^ r dr . (3) 
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The function J(r) contains the information about the structure of the environment ((3 = l/(fcBT))!'! 
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If ReZAuj) f , l3hto r „ . . , 

J[t) = — — / s cotn — — — [cos cot — 1J — z sin lot > dicj , (4) 
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where Z t {uj) = X/{iuC + l/R E )- Note, that at T = one has 

P' (2^) ^ M(2eU) for -► 0, 

P' (2^1) =► ft<5(2eV - 4£ c ) for -> oo. 

This peak structure actually has been seen in experiment!^. The result is said to be correct if the Josephson coupling 
energy Ej — h/(2e)I c is much smaller than E c . By use of the known formula for the critical current I c one gets the 
inequality 

Ec»^A. (5) 

R is the normal tunnel resistance which obeys the relation Rq <C R- Unfortunately, this model contains nearly no 
information about the superconducting electrodes. In contrast to this the quasiparticle current (T > 0) 

oo 

f f eV \ 1 - er$ eV 

(I) qp (V) = J ImJ,( W )P [ T -») T^KJ^ (6) 

— oo 

is expressed in terms of the quasiparticle current amplitude Im/ g which depends in a characteristic way oruoj. In case 
of a symmetric junction the current amplitude lml q (uj) reads at T — according to standard BCS theoryGJ 

!0 for < r -f < 2 

(v/l-(M) 2 ) - ^l-m 2 ) for % > 2 ' 

The symbols E and X stand for the complete elliptic integrals of the first kind0. The definition of P(u>) differs from 
that of P'(uj) (Eq. (0)) by the lack of the factor 4 in front of the function J(r). The dependence of the supercurrent 
on the factor 1^ in Eq. ^ indicates that the supercurrent has something to do with the squared Josephson current 
amplitude Re/ p . This is motivation to express the supercurrent (I) s by means of a perturbation theory of higher 
order in the elementary tunneling Hamiltonian 

H T = H + + H- , H-=Hl, 

H+ = Y, Tir^ct^e** , (8) 

where c; i(T and c r , a stand for quasiparticle annihilation operators of the left and right electrode satisfying anticommu- 
tation relations. The spin is labeled by the subscript a. 

In this way the special features of Cooper pair tunneling (transfer of charges 2e, energy transfer only to the 
environment) arise automatically. In other words we do not consider tunneling particles with charge 2e from the 
beginning. Rather than we start with elementary particles (electrons) with charge le and the supercurrent arises 
as an effect of higher order. Furthermore, the dependence on Re/ p describes the transition to a new branch if the 
supplied energy is able to break Cooper pairs into quasiparticles (eV = 2 A). 

Using the expression for the mean current in nonlinear response theory those parts which correspond to the Cooper 
pair tunneling can be identified and calculated in a systematic way. 

II. RESPONSE THEORY 

The dynamics of a physical system modeled by H = H Q + will be described by the statistical operator p 
satisfying the von Neumann equation which is in the interaction representation (superscript (/)) equivalent to the 
integral equation 
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P (t)^= Po -j / [H^(t'), P ^(t')]dt' . (9) 



i? is the unperturbed part of the Hamiltonian whereas the interaction part Ht reads in the interaction representation 

It is assumed that the interaction is switched on at t = — oo adiabatically. The operator p Q is given by the canonical 
expression 



Po = 



The solution can be found by successive approximation. Now the mean value of the current operator reads 
t 



(J> = 4 / dt'([I^(t),H^(t')]) 



t' t" 



+ / dt ' I dt ' I d ^'([[[ /(/) W'4 7) (0],4 /) (n],4 /) r)])o + ... (io) 



-oo — oo 



In this equation the mean values have to be calculated with respect to p . A term of zeroth order is missing because 
in case of no interaction (tunneling) there is also no current. The first term corresponds to linear response theory and 
leads either to the known quasiparticle tunneling (Eq. (|^)) or to the Josephson current. The second term describes 
the first nonlinear corrections. 



III. THE MODEL 

Now let us apply this theory to tunneling through a SIS junction with environment (Fig. Q). The total Hamiltonian 
reads 

H = H + H T =QV + H res + H T , (11) 

where the tunnel Hamiltonian is given by Eq. (^|). In case ot-superconducting electrodes one can assume that the 
macroscopic phase is already contained in the phase operator $113. Owing to this phase operator tunneling is connected 
with excitations in the electromagnetic environment. The operator H + means e.g. tunneling from left to right in 
contrast to the Hermitian conjugate which describes the reverse process. Ti r are the tunneling matrix elements. The 
basic algebra ruling this approach is the following relationEj 

H± ■ F(Q) = F(Q ±e)-H±, (12) 

where F is an arbitrary function of the junction charge Q. This algebra corresponds to the elementary commutation 
relation 

[Q^)=ie. (13) 

The convention is chosen in such a way that a positive voltage favors tunneling from left to right which reduces the 
junction charge Q by e. 

The reservoir Hamiltonian H res consis±s_of terms corresponding to the left and right electrodes and the environment 
which can be described in standard wayEffl. 

For the calculation of the stationary mean current in terms of Eq. ( |l0| ) a current operator has to be defined. This 
is done in the following way 

I = ~Q = ~[Q, H] - ±(H + HJ) . (14) 
d£ ih in 

Now the quasiparticle current which is the first order term in Eq. ( fl0| ) reads 
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(I) qp = -^Re J dt'([H^(t),H^(t')]) . (15) 



The Cooper pair current is contained in the following terms of second nonvanishing order 

t t' t" 

(J) s = ^Re / dt> [ dt" [ dt"'{([[[H^\t),Hi I \t%H^\t'%H^\n}) 



+ [[[H^(t), (*')], H$> (*")], H^\t"')}) 

+ [[[H^\t),H^\t%H^\t")],H^\t'")]) } . (16) 

An explanation should be given why only correlations with vanishing signature (+ H and their permutations) 

are taken into account. The reason is the following. The separation of the voltage dependence in the correlation 
functions ((. . .} ) without signature zero (e.g. (H + H + H + H + ) ) leads to expressions containing time dependent (t) 
terms. Furthermore, a phase <p remains indeterminated additionally. This corresponds to common Josephson physics 
where one has contributions proportional to sin(2eVt / h + <j) ) and cos(2eVt / h + o )t3. However, in Josephson physics 
the phase <j) is defined by current biasing. Therefore, in our case of voltage-biasing one has to average with respect 
to this phase factor cj> which yields zero. 

By splitting-off the voltage dependence by means of Eq. (|l2|) and using new time variables 



* - t': t' =t' - t"; 



Eq. (Ill) reads 



oo oo oc 

(7) s = pRe JdrJ dr' J dr" [e^ v ^' ^^{r, r>, r") 



+e -^(^") K2 ( T y ;T ") 



IV. CALCULATION OF CORRELATION FUNCTIONS 



Now one has to deal with the three correlation functions K\, «2 and K3. The operators H±(t) in the interaction 
representation can be written as H±(t) exp(±i#(t)) where the operators H±(t) carry only the time dependence with 
respect to the electrodes and the phase dependend operators exp(±i<P(t)) carry those with respect to the environment. 
For instance the function k\ reads therefore 

Ki = (H+(t)H+(t - r)H-(t - t - r')H-(t -t-t' - t")) 

x ^ e »(t) e »(t-r) e -»(t-T-T') e -i*(t-r-T'-r")) o + 7 f urt hcr terms . (18) 

These other terms arise due to the resolution of the interlaced commutators. The decisive step of the identification 
of the contributions which describe Cooper pair tunneling is to reduce the 4-correlators with respect to H into 
2-correlators containing only operators with the same signature, namely 

( J ff + (t 1 ) J ff + (t 2 ) J ff_(t 3 )^-(i4))o = (H+(h)H + (t 2 )) {H_{h)H_{t A )) . 

This decomposition guarantees that only condensate states corresponding to Cooper pairs are taken into account. One 
can prove this from the point of view of the elementary operators c] and c ] . Then the decomposition is equivalent 
to 

(4(t 1 )4(t3)} (c I .^ 2 )c r (t 4 )) (Q(ti)cKi3))o(4^2)c ; t (t 4 )}o 

and one can see that on both banks of the junction only the condensate properties contribute. The terms which have 
been neglected in this decomposition belong to quasiparticle tunneling of higher order and processes including both 
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quasiparticles and Copper pairs. The correlators (H+(ti)H+(t2)) can be expressed by the current amplitude lml p 
in the following waytj 



«±(r) = (H±(t)H±(t - t)) = (H±(t)H ± (0)) o 



9 w 

P C 1 

- — / dwlm/„(w)e-^ r 

2ne J pK ' l-e- h < 



w/(fc fc T) 



(19) 



Note the symmetry k+(t) = k_(t). The phase correlations can also be calculated, for instance by generalizing the 
method presented in Ref. ||. 

Now the correlation function ki(t, t , r") can be written as 

K X {t,t',t") = K + (r) K _(r")e J ( r+4T ' +T ") - k + (-t)k_ (_ t ")^(---4t'-t") 

- K+ (-r) K _(r")e J(3r+4T ' +r " ) + K+ (t)k_ (_ r ") e ^(-3--4r'-r") 
-K + (r)K_(r")e J ^ T+T ") + K + (-r)/e_(-<r")e J < T - T "> 

+ K+ (-r) K _(r")e J(T+T " ) - K+ (r) K _(-r")e J( - T - T " ) . (20) 

For the environment function J(t) see Eq. (0). We remind of the property J(r) = —ilo c t at T = and for Re/Rq — » 
oo. 



V. SUPERCURRENT 

Now the contribution to (J) s coming of «i ((J)a ) can be calculated. Using the definition 

Im/ P (w) t^o 



/(") 



ImJ p (a;)0(cj) , 



where is the unit step function as well as the definition of the function <5 + 



(21) 



S+(x) 



1 

2^ 



z lkx Ak 



the result can be written as 



oo oo 

(I) ( s Kl) = Y Re J duj J dw7H/(w') I 5+(-v-u-u c )S + (-v-u'-u e )8 + (-v-2u c ) 



— oo -co 



-s + 

-8+ 
+ 8+ 
-8+ 
+ 8+ 
+ 8+ 
-8, 



-v + to + lo c )8+ (— v + u) + uj c )d+(—v + 2oj c ) 
-v + w — 3lu c )5 + (— v — to' — lu c )S+(—v — 2lu c ) 
-v — ui + 3lu c )5 + (— v + uj + lu c )5 + (—v + 2uj c ) 

-V — LU + LO c )8+( — V — U) — LJ c )S + ( — v) 

-v + 10 — U! c )S + (—V + to' + lo c )S + (— v) 

-V + LO — LU c )S + ( — V — LU 1 — LO c )8 + ( — v) 

-v — to + lu c )5 + (—v + to' + lo c )5 + (— u)| . 



(22) 



Here, the variable v = eV/Ti is used. The same procedure has to be employed with respect to the terms including the 
other correlation functions «2(t, t',t") and «3(r, t',t"). The real part of the sum of these terms can be calculated 
by means of the Dirac formula 



8 + { x )= l -Ux)+V l - 



(23) 



Finally, at least one integration (e.g. with respect to to', see Eq. (p2|)) can be carried out and one gets after a lenghty 
calculation 
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(I)s(v) 



2e 



OO 

1 r a.^l 



TT UJ — UJ r 



f(v - UJ C ) 2oj c 



f 



+ f(u c ) } 6(v ~ 2uj c ) 

f(w)(u + UJ C 



2ire v(v — 2uj c ) J (u> + v — lu c )(lu — v + u> c ) 



The dash in the integral sign means that one has to take the principal value of the integral. Eq. (| 
result. 



(24) 



is our mam 



VI. DISCUSSION 



Concerning the structure of the current let us concentrate on the case of a high-resistance environment (u> c ^ 0) . 
The opposite case of a low-resistance environment corresponds to the condition oj c = 0. One can make the following 
remarks: 

(i) The current is an antisymmetric function of the applied voltage which reflects the expectation that a reversed 
voltage leads to a reversed current. 

(ii) The current shows a <5-like singularity at 2eV = 4E C corresponding to the fact that the energy 2eV connected 
with the tunneling of a Cooper pair has to be transferred to the environment. Because Cooper pairs live in the 
condensate they cannot absorb this energy. Of course, this singular expression will be smoothed both due to 
finite temperatures and finite environment resistances. 

(iii) There is an additional current contribution which is proportional to f(v — u> c )—. Because of Eq. ( |2l| ) and of the 
known structure of lml p (uj) in standard BCS theory (ImJ p (w) = — ImJ p (— w))til 



ImL(w) = —I c 

7T 



TiUJ 



l-(ff) 



for < < 2 



for ^ > 2 



(25) 



this current contribution only exists if v — uj c < 2A/U or eV < 2 A + E c which is just the condition for the 
onset of the quasiparticle current. Since our approach is based on higher order perturbation theory we are only 
interested in effects which occur in the gap region of quasiparticle tunneling. Therefore, only the first term of 
Eq. ( pi|) has to be considered. 

The integral in the first term of Eq. (Q) reminds of the definition of Re/ p (a; c ) according to the Kramers-Kronig 
relation. The only difference is the O-function in the integrand. Nevertheless, it is reasonable to discuss the case 
uj c < 2A/U corresponding to the subgap region between and the position of the Riedel peak. This means f{oJ c ) = 
and the supercurrent reads for < v < 2A/U + uj c 



7T 

Ye 



(I)s(v 

Fig. ||, where we have plotted the expressions 

OO 

1 f , , JmlJuj') 2 

TT J UJ — UJ TT 



S(v ~ 2uj c ) . 



K (§|) for < ^ < 2 



(26) 



2A 



~A > Z 



(27) 



(see Ref. [□]) and 
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If 



do; 



, Iml p {uj')<d(uj') 



(28) 



G 



respectively, shows for < u < 2A/h that sufficiently far from the position of the Riedel peak the approximation 

oo 

_ I I fa J^L „ lRel p (co c ) (29) 

7T J U! — U! c I 

— oo 

holds which becomes exact for lo c — > 0. Hence, for u c < 2A/% one can write 

' 7 d „^U^,. (30, 

Using this approximation, the supercurrent reads 

(I)s(v) = ^I 2 c {5{v - 2u ) - S(-v - 2u c )} (31) 

8e 

which corresponds for T — > and Re/Rq — > oo exactly to the result of Eq. (j^). It has been shown that this 
formula is valid for < eV < 2A + E c and E c <C 2A There is no contradiction to the inequality (||) because 
2A ^> E c ^> Rq/(8R)A is satisfied provided that the relation Rq <C R holds. However, this condition is just necessary 
for single-charge tunneling because it guarantees that quantum fluctuations can be neglected. Roughly speaking (cf. 



(|29D), formula ( |26[ ) shows indeed that the Josephson current amplitude Kel p (uj c ) determines the strength of the <5-like 
current peak at the eV = 2E C . Note, that this strength becomes singular if the (5-singularity tends to the onset 
position of quasiparticle tunneling because from v = 2lo c — > 2A/H + lo c follows uj c — > 2 A/%. Eq. ( |2^ ) shows that for 
T > there are also current contributions depending on the dissipative part of the Josephson current Im/ p which 
describes pair transfer processes via thermally excited quasiparticles. 

In case of a finite environment resistance the substitution J(t) = —iuj c t in Eq. ( pfj| ) is not possible. The investigation 
has shown that the origin of the resulting (5-function in Eq. ( pl| ) is the integration over t' in Eq. ([l7]). Therefore, 
this integration would indeed generate the function P'(2eV) known from Eq. (||). But there are also functions J(t) 
and J(j") which are modifying the other integrations over r and t" . To sum up it can be said that in this stricter 
approach the dependence on the environment is much more complicated than in the model (0) which leads to Eq. (|^). 

VII. RECONSTRUCTION OF AN EFFECTIVE HAMILTONIAN 

Using some simple assumptions we are going to reconstruct an effective Hamiltonian which leads in first order 
perturbation theory to the same result @. It can also be written as 

(32) 

It will turn out that this effective Hamiltonian corresponds just to the Hamiltonian (|l|). 

Our starting point is the Hamiltonian ( [TT] ) and the assumption that the perturbation term H? can be written as 

Ht = H + + H- = He iw + He- 1 * . 

This ansatz with real constants H seems to be very likely because the factor in front of the (5-function in the super- 
current is also a constant. W is a phase operator which is assumed to obey the commutation relation 

[Q, V] - ike , (33) 

where the constant k is for the time being arbitrary. One gets the algebra 

H±F{Q) = F(Q ± ke)H ± (34) 

and finds using Eq. (|lj) by linear response theory the mean current 




</> s = -|£Re J d f <[fff(t),ffi'>(f)]) 



(35) 
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with 



The time dependence of the phase operator governed by the environment Hamiltonian can be calculated in a standard 
way. Finally, Eq. (pq) reads 



2keH 2 ^ 



(I) s = -~— Re / dre-« feeV/r e fc J(T > - e fc J( - T) , (36) 

h J l j 

o 

where we already know the function J from Eq. (0). By comparing Eq. (^) with Eq. ( |32"|) in the limit T — the 
unknown constants and k can be determined 

fc = 2, H =^-I c = ^- (37) 

which reproduce expression (|l|) with = 2^. The value k = 2 shows that the effective Hamiltonian describes tunneling 
of electron pairs (Cooper pairs). In this way the transition to the effective model corresponds to the transition from 
[Q,$] = ie to [Q,&] = «2e. However, the effective model does not contain dissipation because it depends instead of 
I p (u>) on I c only. 
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FIG. 1. Scheme of the circuit 

FIG. 2. Plot of the expressions according to Eq. (^) (top) and ( pq ) (bottom) in units of I c versus lu in units of A/h in the 
subgap region (0 < u < 2A/h) 
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